home *** CD-ROM | disk | FTP | other *** search
/ Linux Cubed Series 3: Developer Tools / Linux Cubed Series 3 - Developer Tools.iso / devel / lang / lisp / guile-ii.src / guile-ii / guile-src / slib / array.scm < prev    next >
Encoding:
Text File  |  1993-11-02  |  8.1 KB  |  280 lines
  1. ;;;;"array.scm" Arrays for Scheme
  2. ; Copyright (C) 1993 Alan Bawden
  3. ;
  4. ; Permission to copy this software, to redistribute it, and to use it
  5. ; for any purpose is granted, subject to the following restrictions and
  6. ; understandings.
  7. ;
  8. ; 1.  Any copy made of this software must include this copyright notice
  9. ; in full.
  10. ;
  11. ; 2.  Users of this software agree to make their best efforts (a) to
  12. ; return to me any improvements or extensions that they make, so that
  13. ; these may be included in future releases; and (b) to inform me of
  14. ; noteworthy uses of this software.
  15. ;
  16. ; 3.  I have made no warrantee or representation that the operation of
  17. ; this software will be error-free, and I am under no obligation to
  18. ; provide any services, by way of maintenance, update, or otherwise.
  19. ;
  20. ; 4.  In conjunction with products arising from the use of this material,
  21. ; there shall be no use of my name in any advertising, promotional, or
  22. ; sales literature without prior written consent in each case.
  23. ;
  24. ; Alan Bawden
  25. ; MIT Room NE43-510
  26. ; 545 Tech. Sq.
  27. ; Cambridge, MA 02139
  28. ; Alan@LCS.MIT.EDU
  29.  
  30. (require 'record)
  31.  
  32. ;(declare (usual-integrations))
  33.  
  34. (define array:rtd
  35.   (make-record-type "Array"
  36.     '(indexer        ; Must be a -linear- function!
  37.       shape        ; Inclusive bounds: ((lower upper) ...)
  38.       vector        ; The actual contents
  39.       )))
  40.  
  41. (define array:indexer (record-accessor array:rtd 'indexer))
  42. (define array-shape (record-accessor array:rtd 'shape))
  43. (define array:vector (record-accessor array:rtd 'vector))
  44.  
  45. (define array? (record-predicate array:rtd))
  46.  
  47. (define (array-rank obj)
  48.   (if (array? obj) (length (array-shape obj)) 0))
  49.  
  50. (define (array-dimensions ra)
  51.   (map (lambda (ind) (if (zero? (car ind)) (cadr ind) ind))
  52.        (array-shape ra)))
  53.  
  54. (define array:construct
  55.   (record-constructor array:rtd '(shape vector indexer)))
  56.  
  57. (define (array:compute-shape specs)
  58.   (map (lambda (spec)
  59.      (cond ((and (integer? spec)
  60.              (< 0 spec))
  61.         (list 0 (- spec 1)))
  62.            ((and (pair? spec)
  63.              (pair? (cdr spec))
  64.              (null? (cddr spec))
  65.              (integer? (car spec))
  66.              (integer? (cadr spec))
  67.              (<= (car spec) (cadr spec)))
  68.         spec)
  69.            (else (slib:error "array: Bad array dimension: " spec))))
  70.        specs))
  71.  
  72. (define (make-array initial-value . specs)
  73.   (let ((shape (array:compute-shape specs)))
  74.     (let loop ((size 1)
  75.            (indexer (lambda () 0))
  76.            (l (reverse shape)))
  77.       (if (null? l)
  78.       (array:construct shape
  79.                (make-vector size initial-value)
  80.                (array:optimize-linear-function indexer shape))
  81.       (loop (* size (+ 1 (- (cadar l) (caar l))))
  82.         (lambda (first-index . rest-of-indices)
  83.           (+ (* size (- first-index (caar l)))
  84.              (apply indexer rest-of-indices)))
  85.         (cdr l))))))
  86.  
  87. (define (make-shared-array array mapping . specs)
  88.   (let ((new-shape (array:compute-shape specs))
  89.     (old-indexer (array:indexer array)))
  90.     (let check ((indices '())
  91.         (bounds (reverse new-shape)))
  92.       (cond ((null? bounds)
  93.          (array:check-bounds array (apply mapping indices)))
  94.         (else
  95.          (check (cons (caar bounds) indices) (cdr bounds))
  96.          (check (cons (cadar bounds) indices) (cdr bounds)))))
  97.     (array:construct new-shape
  98.              (array:vector array)
  99.              (array:optimize-linear-function
  100.                (lambda indices
  101.              (apply old-indexer (apply mapping indices)))
  102.                new-shape))))
  103.  
  104. (define (array:in-bounds? array indices)
  105.   (let loop ((indices indices)
  106.          (shape (array-shape array)))
  107.     (if (null? indices)
  108.     (null? shape)
  109.     (let ((index (car indices)))
  110.       (and (not (null? shape))
  111.            (integer? index)
  112.            (<= (caar shape) index (cadar shape))
  113.            (loop (cdr indices) (cdr shape)))))))
  114.  
  115. (define (array:check-bounds array indices)
  116.   (or (array:in-bounds? array indices)
  117.       (slib:error "array: Bad indices for " array indices)))
  118.  
  119. (define (array-ref array . indices)
  120.   (array:check-bounds array indices)
  121.   (vector-ref (array:vector array)
  122.           (apply (array:indexer array) indices)))
  123.  
  124. (define (array-set! array new-value . indices)
  125.   (array:check-bounds array indices)
  126.   (vector-set! (array:vector array)
  127.            (apply (array:indexer array) indices)
  128.            new-value))
  129.  
  130. (define (array-in-bounds? array . indices)
  131.   (array:in-bounds? array indices))
  132.  
  133. ; Fast versions of ARRAY-REF and ARRAY-SET! that do no error checking,
  134. ; and don't cons intermediate lists of indices:
  135.  
  136. (define (array-1d-ref a i0)
  137.   (vector-ref (array:vector a) ((array:indexer a) i0)))
  138.  
  139. (define (array-2d-ref a i0 i1)
  140.   (vector-ref (array:vector a) ((array:indexer a) i0 i1)))
  141.  
  142. (define (array-3d-ref a i0 i1 i2)
  143.   (vector-ref (array:vector a) ((array:indexer a) i0 i1 i2)))
  144.  
  145. (define (array-1d-set! a v i0)
  146.   (vector-set! (array:vector a) ((array:indexer a) i0) v))
  147.  
  148. (define (array-2d-set! a v i0 i1)
  149.   (vector-set! (array:vector a) ((array:indexer a) i0 i1) v))
  150.  
  151. (define (array-3d-set! a v i0 i1 i2)
  152.   (vector-set! (array:vector a) ((array:indexer a) i0 i1 i2) v))
  153.  
  154. ; STOP!  Do not read beyond this point on your first reading of
  155. ; this code -- you should simply assume that the rest of this file
  156. ; contains only the following single definition:
  157. ;
  158. ;   (define (array:optimize-linear-function f l) f)
  159. ;
  160. ; Of course everything would be pretty inefficient if this were really the
  161. ; case, but it isn't.  The following code takes advantage of the fact that
  162. ; you can learn everything there is to know from a linear function by
  163. ; simply probing around in its domain and observing its values -- then a
  164. ; more efficient equivalent can be constructed.
  165.  
  166. (define (array:optimize-linear-function f l)
  167.   (let ((d (length l)))
  168.     (cond
  169.      ((= d 0)
  170.       (array:0d-c (f)))
  171.      ((= d 1)
  172.       (let ((c (f 0)))
  173.     (array:1d-c0 c (- (f 1) c))))
  174.      ((= d 2)
  175.       (let ((c (f 0 0)))
  176.     (array:2d-c01 c (- (f 1 0) c) (- (f 0 1) c))))
  177.      ((= d 3)
  178.       (let ((c (f 0 0 0)))
  179.     (array:3d-c012 c (- (f 1 0 0) c) (- (f 0 1 0) c) (- (f 0 0 1) c))))
  180.      (else
  181.       (let* ((v (map (lambda (x) 0) l))
  182.          (c (apply f v)))
  183.     (let loop ((p v)
  184.            (old-val c)
  185.            (coefs '()))
  186.       (cond ((null? p)
  187.          (array:Nd-c* c (reverse coefs)))
  188.         (else
  189.          (set-car! p 1)
  190.          (let ((new-val (apply f v)))
  191.            (loop (cdr p)
  192.              new-val
  193.              (cons (- new-val old-val) coefs)))))))))))
  194.  
  195. ; 0D cases:
  196.  
  197. (define (array:0d-c c)
  198.   (lambda () c))
  199.  
  200. ; 1D cases:
  201.  
  202. (define (array:1d-c c)
  203.   (lambda (i0) (+ c i0)))
  204.  
  205. (define (array:1d-0 n0)
  206.   (cond ((= 1 n0) +)
  207.     (else (lambda (i0) (* n0 i0)))))
  208.  
  209. (define (array:1d-c0 c n0)
  210.   (cond ((= 0 c) (array:1d-0 n0))
  211.     ((= 1 n0) (array:1d-c c))
  212.     (else (lambda (i0) (+ c (* n0 i0))))))
  213.  
  214. ; 2D cases:
  215.  
  216. (define (array:2d-0 n0)
  217.   (lambda (i0 i1) (+ (* n0 i0) i1)))
  218.  
  219. (define (array:2d-1 n1)
  220.   (lambda (i0 i1) (+ i0 (* n1 i1))))
  221.  
  222. (define (array:2d-c0 c n0)
  223.   (lambda (i0 i1) (+ c (* n0 i0) i1)))
  224.  
  225. (define (array:2d-c1 c n1)
  226.   (lambda (i0 i1) (+ c i0 (* n1 i1))))
  227.  
  228. (define (array:2d-01 n0 n1)
  229.   (cond ((= 1 n0) (array:2d-1 n1))
  230.     ((= 1 n1) (array:2d-0 n0))
  231.     (else (lambda (i0 i1) (+ (* n0 i0) (* n1 i1))))))
  232.  
  233. (define (array:2d-c01 c n0 n1)
  234.   (cond ((= 0 c) (array:2d-01 n0 n1))
  235.     ((= 1 n0) (array:2d-c1 c n1))
  236.     ((= 1 n1) (array:2d-c0 c n0))
  237.     (else (lambda (i0 i1) (+ c (* n0 i0) (* n1 i1))))))
  238.  
  239. ; 3D cases:
  240.  
  241. (define (array:3d-01 n0 n1)
  242.   (lambda (i0 i1 i2) (+ (* n0 i0) (* n1 i1) i2)))
  243.  
  244. (define (array:3d-02 n0 n2)
  245.   (lambda (i0 i1 i2) (+ (* n0 i0) i1 (* n2 i2))))
  246.  
  247. (define (array:3d-12 n1 n2)
  248.   (lambda (i0 i1 i2) (+ i0 (* n1 i1) (* n2 i2))))
  249.  
  250. (define (array:3d-c12 c n1 n2)
  251.   (lambda (i0 i1 i2) (+ c i0 (* n1 i1) (* n2 i2))))
  252.  
  253. (define (array:3d-c02 c n0 n2)
  254.   (lambda (i0 i1 i2) (+ c (* n0 i0) i1 (* n2 i2))))
  255.  
  256. (define (array:3d-c01 c n0 n1)
  257.   (lambda (i0 i1 i2) (+ c (* n0 i0) (* n1 i1) i2)))
  258.  
  259. (define (array:3d-012 n0 n1 n2)
  260.   (cond ((= 1 n0) (array:3d-12 n1 n2))
  261.     ((= 1 n1) (array:3d-02 n0 n2))
  262.     ((= 1 n2) (array:3d-01 n0 n1))
  263.     (else (lambda (i0 i1 i2) (+ (* n0 i0) (* n1 i1) (* n2 i2))))))
  264.  
  265. (define (array:3d-c012 c n0 n1 n2)
  266.   (cond ((= 0 c) (array:3d-012 n0 n1 n2))
  267.     ((= 1 n0) (array:3d-c12 c n1 n2))
  268.     ((= 1 n1) (array:3d-c02 c n0 n2))
  269.     ((= 1 n2) (array:3d-c01 c n0 n1))
  270.     (else (lambda (i0 i1 i2) (+ c (* n0 i0) (* n1 i1) (* n2 i2))))))
  271.  
  272. ; ND cases:
  273.  
  274. (define (array:Nd-* coefs)
  275.   (lambda indices (apply + (map * coefs indices))))
  276.  
  277. (define (array:Nd-c* c coefs)
  278.   (cond ((= 0 c) (array:Nd-* coefs))
  279.     (else (lambda indices (apply + c (map * coefs indices))))))
  280.